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Nomenclature 

Young's modulus of skin and stringer flange 
Young's modulus of adhesive 

bending rigidity of skin and stringer flange 
shear modulus of adhesive 
thickness of skin, flange, and adhesive 
flange width 

vertical shear force applied to skin and flange 
moment applied to skin and flange 
inplane force applied to skin and flange 
( t x +t c ) /D x , ( t 2 +t c ) /D 2 

inplane displacement of skin and flange 
vertical deflection of skin and flange 
eigenvector 

eigenvalue (LAM(i) in program output) 

V X ( i) (MU(i) in program output) 
coefficient 

inplane stress in skin and flange (SIG1,SIG2 in 
program output) 

peel stress in adhesive (SIGC in program output) 
shear stress in adhesive (TAUC in program output) 
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1. INTRODUCTION 


Computer program PSTRESS, which carries out an analysis of two 
bonded plates, was developed as part of NASA contract 
NAS1-17970. The goal of the NASA project is the development of 
key technology for the application of composite materials to 
high aspect ratio transport wing structure. The objectives of 
this effort are the design, manufacture, and test of cover panel 
structure for a baseline composite wing that will meet all 
strength, stiffness, aeroelastic, and damage tolerance 
requirements with maximum weight savings and minimum cost. 

Test results have shown that impact in the area of a T-shaped 
stringer bonded to a composite skin can result in a debond. 
PSTRESS is a computer program developed to investigate the 
influence of various parameters on the adhesive stresses in such 
a bonded composite panel. It is an interactive program written 
in the FORTRAN 77 language. Input includes properties of the 
skin and stringer flange, loading on the skin and stringer 
flange, and properties of the adhesive. Output includes 
functions for the calculation of vertical and inplane 
displacements, inplane stress in the skin and stringer flange, 
and peel and shear stresses in the adhesive layer. These 
quantities are calculated at twenty points across the width of 
the flange. These numbers are presented in the output in 
tabular form. 



2. ANALYTICAL DEVELOPMENT 


Figure 1 shows the model used in the PSTRESS analysis. 

Subscripts 1, 2, and c refer to properties of the skin, stringer 
flange, and adhesive, respectively. A complete description of 
the analysis can be found in Reference 1. Briefly, a system of 
two fourth order and two second order differential equations 
with associated boundary equations are obtained using the 
principle of minimum potential energy. PSTRESS solves for the 
eigenvalues and eigenvectors of the system of equations and 
applies the boundary conditions to solve for the coefficients of 
functions for w^_ , W 2 , u^, and u 2 . Beam bending, linear 
elastic behavior of the material, and symmetry about the plane 
through the blade of the stringer are assumed throughout the 
analysis . 

Peel stress and shear stress in the adhesive are determined from 
the vertical deflections (W]_, W 2 ) and inplane displacements 
(U]_, u 2 ). Equations la and lb show the relations for these 
stresses in terms of u^, u 2 , and w 2 . 

t c =G c ( 2 ( u^~U 2 ) tD]_T^dw 1 /dx+D 2 T 2 dW 2 /dx) /2t c (la) 

a c =E c ( Wi“W 2 ) /t c (lb) 


PSTRESS calculates eigenvalues ( A(i)'s) and eigenvectors 
(U(ij)'s). Boundary conditions are then applied and the necessary 
coefficients of the solutions are determined. Figure 2 shows a 
simplified flow chart of the program. 
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Figure 2 

Flow Chart of PSTRESS 
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The solutions for the system of differential equations take the 
following form: 

Case 1 

Inplane force (F^ or F 2 ) applied. 

w 1 =[C(l)e u(1)x +C(2)e' u(1)x +C(3)e u(2)x +C(4)e -u(2)x +C(5)e u(3)x 
+C(6)e -u(3)x +C(7)e u(4)x +C(8)e _u(4)x +C(9)x+C(10) ]/D x 

w 2 =[U(21){C(l)e u(1)x +C(2)e -u(1)x }+U(22){C(3)e u ' (2)x +C(4)e“ u(2)x } + 

U(23){C(5)e u(3)x +C(6)e" u<3)x >+U(24) {C(7)e u(4)x +C(8)e _u(4)x }]/D 2 
+[C(9)x+C(10) ]/D x 

u 1 =[U(51){C(l)e u(1)x -C(2)e" u(1)x }/u(l)+U(52){C(3)e ,l(2)x - 
C(4)e _u(2)x }/u(2)+U(53){C(5)e u(3)x -C(6)e -u(3)x }/u(3)+ 
U(54){C(7)e u(4)x -C(8)e _u(4)x }/u(4)+C(ll)x+C(13) 1/Eit! 

u 2 =(-U(51){C(l)e u(1)x -C(2)e -ull)x }/u(l)-U(52){C(3)e ,l(2)x - 
C(4)e" u(2)x }/u(2)-U(53){C(5)e u(3)x -C(6)e _w<3)x }/u(3)- 
U(54) {C(7)e u(4)x -C(8)e -ti(4)x }/u(4)+C(12)x+C(14) ]/E 2 t 2 

(2a-d) 

Case 2 

Inplane force not applied (Fi=F2=0 ) . 

w 1 =tC(l)e u(1)x +C(2)e- li(1)x +C(3)e u(2)x + c(4)e-' l(2)x + C(5)e ll(3)x 
+C(6)e~ u(3)x +C(7)x 3 +C(8)x 2 +C(9)x+C(10) ]/D! 

w 2 =[U(21){C(l)e u(1)x +C(2)e" u(1)x }+U(22){C(3)e ll(2)x +C(4)e -u(2)x } + 

U(23){C(5)e u(3)x +C(6)e' u(3)x >]/D 2 +[C(7)x 3 +C(8)x 2 +C(9)x+ 

CdOlJ/D! 
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u 1 =[U(51){C(l)e u(1)x -C(2)e' u(1)x }/u(l)+U(52){C(3)e u(2)x - 
C(4)e -U(2)X }/u(2)+U(53){C(5)e W(3)X -C(6)e~ u(3,X }/u(3)- 
3E 1 t 1 E 2 t 2 (T 1 +T 2 D 2 /D 1 ) C ( 7 ) X 2 / ( E 1 t 1 +E 2 t 2 ) +C ( 11 ) X+C ( 13 ) ] /E^ 

u 2 =[-U(51){C(l)e u(1)x -C(2)e _li(:L)x }/u(l)-U(52){C(3)e u(2)x - 

C(4)e' u(2)X }/u(2)-U(53){C(5)e u(3)x -C(6)e -u(3)x }/u(3)+ 

3E 1 t 1 E 2 t 2 (T 1 +T 2 D 2 /D 1 )C(7)x 2 /(E 1 t 1 +E 2 t 2 )+C(12)x+C( 14) ]/E 2 t 2 

( 3a-d) 
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3. PROGRAM INPUT AND OUTPUT 


Program Input 

The operator is prompted for all input. Necessary input for the 
program is as follows: 


El, 

Dl, 

and 

fc l : 

modulus, bending rigidity, 
the skin 

and 

thickness 

of 

E 2< 

»2' 

and 

t 2 : 

modulus, bending rigidity, 
the stringer flange 

and 

thickness 

of 

E C' 

G c 

and 

fc c : 

extensional modulus, shear 
thickness of the adhesive 

modulus , and 



L : width of stringer flange 

P^, M]_ , and F]_: running vertical shear force, running moment, 

and running inplane force applied to the skin 
P 2 , m 2' an ^ f 2 : running vertical shear force, running moment, 

and running inplane force applied to the stringer 
flange 

Program Output 

All output is displayed to the screen and is also written to a 
file entitled "PSTRESS OUTPUT". The following quantities are 
output : 

Echo of input quantities 

Eigenvalues and their square roots: A(i) ! s and u(i)’s 

Eigenvectors: w^ (U(l,j ) ) , W 2 (U(2,j ) ) , d 2 w^/dx 2 

( U( 3 , j ) ) , d 2 w 2 /dx 2 ( U( 4 , j ) ) , and u x (U(5,j) ) . 
Coefficients: C(i)’s 

Deflections: w^, w 2 , U]_, and u 2 for twenty points 

across flange width 

Stresses: a c , t c , a^, and a 2 for twenty points 

across flange width 
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4. SAMPLE PROBLEMS 


The following pages contain the input and output for analyses 
using the data given below: 

£2=4713 ksi (32.473 GPa) 

D2=3685 lb-in (416 N-m) 
t 2 =. 2232 in (5.7 mm) 

E c =500 ksi (3.445 GPa) 
t c =.005 in (.127 mm) 

Case 1 

P2=100 Ibs/in (17.5 kN/m) P 2 =0 

M2=80 lb-in/in (355.6 N-m/m) M 2 =0 

F;j_=16Q0 lbs/in (280 kN/m) F 2 =0 

Case 2 

P 1 =100 lbs/in (17.5 kN/m) P 2 =0 

M2=80 lb-in/in (355.6 N-m/m) M 2 =0 

F i=0 F 2 =0 

The ouput is plotted for each case and shown in Figures 3 
through 11. For comparison, the results of a NASTRAN finite 
element analysis of an equivalent model are also shown in 
figures 3, 4, 8, and 9. 


E 2 =4713 ksi (32.473 GPa) 
D 2 =357 lb-in (40.3 N-m) 
t 2 =. 1116 in (2.8 mm) 
G c =45 ksi (.310 GPa) 
L=1.2 in (30.48 mm) 
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PEEL STRESS 


INPUT PROPERTIES OF SKIN: El, D1 , T1 

7 

4713000.3685. . 2232 

INPUT PROPERTIES OF FLANGE: E2 , D2 , T2 

7 

4713000.357. . 1116 

INPUT ADHESIVE PROPERTIES: EC, GC , TC 

7 

500000.45000. . 005 

INPUT FLANGE WIDTH: L 

7 

i.2 

INPUT SHEAR, MOMENT, AND FORCE APPLIED TO SKIN 

7 

100,80,1600 

INPUT SHEAR, MOMENT, AND FORCE APPLIED TO FLANGE 

7 

0 , 0,0 


DO YOU WANT TO SEE EIGENVALUES, EIGENVECTORS, AND COEFFICIENTS ( Y/N ) ? 
"NO" IS DEFAULT, 
n 


DEFLECTIONS 


X 

W1 

W2 

U1 

U2 

0.0000 

1 . 367E-05 

1 . 381E-34 

1.050E-35 

8 . 999E-36 

0.0188 

1.508E-05 

2 . 555E-06 

1 . 211E-05 

3 . 281E-05 

0.0375 

1 . 932E-05 

9 . 278E-06 

2 . 426E-05 

6 . 555E-05 

0.0750 

3 . 609E-05 

3 . 147E-05 

4 . 875E-05 

1 . 307E-04 

0.1125 

6 . 347E-05 

6 . 247E-05 

7 . 356E-05 

1 . 951E-04 

0.1500 

1 . 008E-04 

1 . 014E-04 

9 . 875E-05 

2 . 588E-04 

0.2250 

2 . 026E-04 

2 . 033E-04 

1.504E-04 

3.836E-04 

0.3000 

3 . 369E-04 

3.371E-04 

2 . 039E-04 

5.047E-04 

0.4200 

6 . 106E-04 

6 . 106E-04 

2 . 938E-04 

6 . 900E-04 

0.5400 

9 . 457E-04 

9 . 458E-04 

3 . 891E-04 

8 . 644E-04 

0.6600 

1.331E-03 

1 . 331E-03 

4 . 900E-04 

1 . 028E-03 

0.7800 

1.756E-03 

1 . 757E-03 

5 . 967E-04 

1 . 179E-03 

0.9000 

2 . 215E-03 

2 . 216E-03 

7 . 101E-04 

1 . 318E-03 

0.9750 

2.518E-03 

2 . 522E-03 

7 . 855E-04 

1.395E-03 

1.0500 

2 . 839E-03 

2 . 853E-03 

8 . 663E-04 

1 . 462E-03 
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1.0875 

3 . 011E-03 

3 . 027E-03 

9 . 097E-04 

1 . 489E-03 

1.1250 

3 . 194E-03 

3 . 203E-03 

9 . 560E-04 

1 . 510E-03 

1.1625 

3 . 394E-03 

3 . 371E-03 

1 . 006E-03 

1 . 525E-03 

1.1812 

3 . 503E-03 

3 . 450E-03 

1 . 032E-03 

1 . 529E-03 

1.2000 

3 . 618E-03 

3 . 526E-03 

1 . 060E-03 

1 . 530E-03 



STRESSES 



X 

SIGC 

TAUC 

SIGl 

SIG2 

0.0000 

1. 367E+03 

-1 . 454E-28 

3.044E+03 

8. 250E+03 

0.0188 

1.253E+03 

4 . 787E+01 

3 . 048E+03 

8 . 240E+03 

0.0375 

1 . 004E+03 

9 . 478E+01 

3 . 060E+03 

8 . 217E+03 

0.0750 

4 . 615E+02 

1 . 806E+02 

3.096E+03 

8 . 145E+03 

0.1125 

1 . 003E+02 

2 . 497E+02 

3 . 141E+03 

8 . 055E+03 

0.1500 

-5 . 543E+01 

3 . 002E+02 

3 . 191E+03 

7 . 955E+03 

0.2250 

-6 . 465E+01 

3 . 573E+02 

3 . 303E+03 

7 . 730E+03 

0.3000 

-1 . 612E+01 

3 . 815E+02 

3 . 427E+03 

7 . 483E+03 

0.4200 

-7 . 243E-01 

3 . 960E+02 

3 . 636E+03 

7 . 065E+03 

0.5400 

-2.991E+00 

4 . 040E+02 

3 . 851E+03 

6 . 634E+03 

0.6600 

-1 . 139E+01 

4 . 216E+02 

4 . 072E+03 

6 . 192E+03 

0.7800 

-2 . 808E+01 

4 . 791E+02 

4 . 312E+03 

5 . 713E+03 

0.9000 

-7 . 832E+01 

6 . 626E+02 

4 . 611E+03 

5.115E+03 

0.9750 

-4 . 089E+02 

9 . 444E+02 

4 . 883E+03 

4 . 571E+03 

1.0500 

-1 . 345E+03 

1.551E+03 

5 . 311E+03 

3 . 715E+03 

1.0875 

-1 . 650E+03 

2 . 095E+03 

5 . 621E+03 

3 . 094E+03 

1.1250 

-8 . 826E+02 

2 . 903E+03 

6 . 019E+03 

2 . 299E+03 

1.1625 

2 . 340E+03 

4 . 061E+03 

6 . 523E+03 

1 . 292E+03 

1.1812 

5 . 308E+03 

4 . 792E+03 

6 . 824E+03 

6 . 890E+02 

1.2000 

9 . 243E+03 

5 . 629E+03 

7 . 168E+03 

1 . 272E-28 


ANOTHER RUN WITH DIFFERENT LOADS (Y/N)? 

Y 

INPUT SHEAR, MOMENT, AND FORCE APPLIED TO SKIN 
7 

100,80,0 

INPUT SHEAR, MOMENT, AND FORCE APPLIED TO FLANGE 
7 

0 , 0,0 


DO YOU WANT TO SEE EIGENVALUES, EIGENVECTORS, AND COEFFICIENTS { Y/N ) ? 
"NO” IS DEFAULT. 

Y 

EIGENVALUES AND THEIR SQUARE ROOTS: 

LAM{1) 2 . 4945642E+01 +5 . 4952145E+02i 

MU { 1 ) 1.6956312E+01 +1 . 6204038E+01i 

LAM( 2 ) 2.4945642E+01 -5 . 4952145E+02i 

MU( 2 ) 1.6956312E+01 -1 . 6204038E+01i 
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LAM( 3 ) 

9 . 3258122E+01 

+0 

MU ( 3 ) 

9 . 6570245E+00 

+0, 

LAM( 4 ) 

O.OOOOOOOE+OO 

+0 , 

MU ( 4 ) 

O.OOOOOOOE+OO 

+0, 

EIGENVECTORS: 


U(l,l) 

1 . OOOOOOOE+OO 

+0, 

U( 1 , 2 ) 

1 . OOOOOOOE+OO 

+0, 

U( 1 , 3 ) 

1. OOOOOOOE+OO 

+0, 

U( 1 , 4 ) 

O.OOOOOOOE+OO 

+0, 

U( 2 , 1 ) 

“8 . 5764511E-01 

+3, 

U(2,2) 

-8 . 5764511E-01 

-3, 

U( 2 , 3 ) 

1 . 0521675E-01 

+0. 

U( 2 , 4 ) 

O.OOOOOOOE+OO 

+0, 

U( 3 , 1 ) 

2 . 4945642E+01 

+5. 

U( 3 , 2 ) 

2 . 4945642E+01 

-5. 

U( 3 , 3 ) 

9 . 3258122E+01 

+0, 

U( 3 , 4 ) 

O.OOOOOOOE+OO 

+0 , 

U(4,l) 

-1 . 9056303E+02 

-4. 

U( 4 , 2 ) 

-1 . 9056303E+02 

+4 , 

U( 4 , 3 ) 

9 . 8123162E+00 

+0 . 

U ( 4,4) 

O.OOOOOOOE+OO 

+0. 

U( 5 , 1 ) 

-9 . 6065846E+02 

+4, 

U( 5 , 2 ) 

-9 . 6065846E+02 

-4. 

U( 5 , 3 ) 

5 . 9785653E+02 

+0. 

U( 5 , 4 ) 

O.OOOOOOOE+OO 

+0, 

COEFFICIENTS ARE: 


C( 1) 

2 . 0393568E-11 

-1. 

C( 2) 

1 . 7996557E-03 

+3. 

C( 3) 

2 . 0393568E-11 

+1 . 

C( 4) 

1 . 7996557E-03 

-3. 

C( 5) 

5 . 6962153E-06 

-4 . 

C( 6) 

7 . 9183185E-02 

+2 . 

C( 7) 

-4 . 2461834E+00 

“8 . 

C( 8) 

2 . 5477184E+01 

+3. 

C( 9) 

7 . 1521850E-01 

+1. 

C( 10 ) 

-3 . 2484452E-02 

-1. 

C(ll) 

-8 . 3589469E+02 

-1. 

C ( 1 2 ) 

8 . 3589469E+02 

+1. 

C ( 13 ) 

4 . 5505142E+00 

+1 . 

C ( 14 ) 

-4 . 5505142E+00 

-1. 


OOOOOOOE+OOi 

OOOOOOOE+OOi 

OOOOOOOE+OOi 

OOOOOOOE+OOi 


OOOOOOOE+OOi 

OOOOOOOE+OOi 

OOOOOOOE+OOi 

OOOOOOOE+OOi 

0784698E-01i 

0784698E-0H 

OOOOOOOE+OOi 

OOOOOOOE+OOi 

4952145E+02i 

4952145E+02i 

OOOOOOOE+OOi 

OOOOOOOE+OOi 

6361495E+02i 

6361495E+02i 

OOOOOOOE+OOi 

OOOOOOOE+OOi 

9828957E+02i 

9828957E+02i 

OOOOOOOE+OOi 

OOOOOOOE+OOi 


6118227E-1H 

4075299E-03i 

6118227E-1H 

4075299E-03i 

3274595E-4H 

1439069E-34i 

9683208E-33i 

2920959E-32i 

3226020E-33i 

4787472E-34i 

0801215E-30i 

0801215E-30i 

2695253E-32i 

2695253E-32i 
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DEFLECTIONS 


X 

W1 

W2 

U1 

U2 

0.0000 

1 . 365E-05 

-2.828E-36 

1 . 494E-37 

-4 . 687E-38 

0.0188 

1 . 626E-05 

3 . 748E-06 

-1 . 417E-05 

2 . 834E-05 

0.0375 

2 . 407E-05 

1.405E-05 

-2 . 831E-05 

5 . 662E-05 

0.0750 

5 . 514E-05 

5 . 056E-05 

-5 . 639E-05 

1 . 128E-04 

0.1125 

1 . 064E-04 

1 . 054E-04 

-8 . 415E-05 

1 . 683E-04 

0.1500 

1 . 771E-04 

1 . 776E-04 

-1 . 115E-04 

2 . 230E-04 

0.2250 

3 . 742E-04 

3 . 748E-04 

-1 . 649E-04 

3 . 299E-04 

0.3000 

6 . 417E-04 

6 . 419E-04 

-2 . 164E-04 

4 . 329E-04 

0.4200 

1 . 207E-03 

1 . 207E-03 

-2 . 944E-04 

5 . 889E-04 

0.5400 

1 . 931E-03 

1 . 931E-03 

-3 . 668E-04 

7 . 336E-04 

0.6600 

2.801E-03 

2 . 801E-03 

-4 . 334E-04 

8 . 668E-04 

0.7800 

3 . 805E-03 

3 . 805E-03 

-4 . 940E-04 

9 . 880E-04 

0.9000 

4 . 935E-03 

4 . 936E-03 

-5 . 479E-04 

1 . 096E-03 

0.9750 

5 . 704E-03 

5 . 707E-03 

-5 . 775E-04 

1 . 155E-03 

1.0500 

6 . 522E-03 

6 . 533E-03 

-6 . 024E-04 

1 . 205E-03 

1.0875 

6 . 953E-03 

6 . 966E-03 

-6 . 125E-04 

1 . 225E-03 

1.1250 

7 . 401E-03 

7 . 408E-03 

-6 . 203E-04 

1 . 241E-03 

1.1625 

7 . 870E-03 

7 . 851E-03 

-6 . 255E-04 

1 . 251E-03 

1.1812 

8 . 114E-03 

8 . 070E-03 

-6 . 269E-04 

1 . 254E-03 

1.2000 

8 . 365E-03 

8 . 287E-03 

-6 . 274E-04 

1 . 255E-03 



STRESSES 



X 

SIGC 

TAUC 

SIG1 

SIG2 

0.0000 

1.365E+03 

7 . 889E-31 

-3 . 564E+03 

7 . 127E+03 

0.0188 

1 . 251E+03 

4 . 881E+01 

-3 . 559E+03 

7 . 118E+03 

0.0375 

1 . 002E+03 

9 . 664E+01 

-3 . 547E+03 

7 . 094E+03 

0.0750 

4 . 589E+02 

1 . 843E+02 

-3 . 510E+03 

7 . 021E+03 

0.1125 

9 . 747E+01 

2 . 551E+02 

-3 . 465E+03 

6 . 929E+03 

0.1500 

-5 . 827E+01 

3 . 073E+02 

-3 . 414E+03 

6 . 827E+03 

0.2250 

-6 . 717E+01 

3 . 675E+02 

-3 . 299E+03 

6 . 597E+03 

0.3000 

-1 . 826E+01 

3 . 943E+02 

-3 . 171E+03 

6 . 342E+Q3 

0.4200 

-2.261E+00 

4 . 120E+02 

-2 . 954E+03 

5 . 908E+03 

0.5400 

-3 . 469E+00 

4 . 214E+02 

-2 . 730E+03 

5 . 460E+03 

0.6600 

-9. 155E+00 

4 . 361E+02 

-2 . 500E+03 

5 . 000E+03 

0.7800 

-1. 996E+01 

4 . 782E+02 

-2 . 256E+03 

4 . 512E+03 

0.9000 

-4 . 764E+01 

6 . 082E+02 

-1 . 970E+03 

3 . 939E+03 

0.9750 

-3 . 025E+02 

8 . 056E+02 

-1 . 730E+03 

3 . 460E+03 

1.0500 

-1 . 084E+03 

1 . 232E+03 

-1 . 379E+03 

2 . 758E+03 

1.0875 

-1 . 371E+03 

1 . 619E+03 

-1 . 135E+03 

2 . 270E+03 

1.1250 

-7 . 798E+02 

2 . 197E+03 

-8 . 315E+02 

1 . 663E+03 

1.1625 

1 . 901E4-03 

3 . 030E+03 

-4 . 593E+02 

9 . 187E+02 

1.1812 

4 . 435E+03 

3 . 556E+03 

-2 . 428E+02 

4 . 856E+02 

1.2000 

7 . 860E+03 

4 . 153E+03 

2 . 631E-03 

-5 . 261E-03 
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ANOTHER RUN WITH DIFFERENT LOADS (Y/N)? 
n 

SEND RESULTS TO PRINTER ( Y/N ) ? 
n 
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Figure 3. Peel stress in adhesive. Case 
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Figure 4. Peel stress in adhesive. Case 2. 
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Figure 5. Shear stress in adhesive. 
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Figure 6. Inplane stress in skin laminate. 
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Figure 7. Inplane stress in flange laminate. 



Figure 8. Out-of-plane displacement of skin laminate. 
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Figure 9. Out-of-plane displacement of flange laminate. 
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Figure 10. Inplane displacement of skin laminate. 




5. SUMMARY 


PSTRESS is an interactive computer program which calculates the 
deflections and stresses of a bonded skin/stringer combination 
with applied moments, inplane forces, and vertical shear 
forces. The results are within a few percent of those obtained 
by NASTRAN finite element analyses of equivalent structures. 
Using PSTRESS it is possible to evaluate the effects of loading 
and material parameters on the adhesive peel and shear stresses. 


- 23 - 


6. REFERENCES 


1. Barkey, D. A., Madan, R. C., and Sutton, J. 0. "Analytical 
Approach to Peel Stresses in Bonded Composite Stiffened 
Panels", 20th Midwestern Mechanics Conference on Composites, 
Purdue University, 1987, (also available as Douglas Paper 
No. 7907) . 


- 24 - 


APPENDIX 
PSTRESS Listing 


- 25 - 


**** P S T R E S S **** 


C 

c 

C **** PROGRAM TO COMPUTE RESPONSE OF BONDED **** 

C **** STRINGER/ SKIN COMBINATION **** 

C 

C **** ANALYSIS AND PROGRAM BY DEREK A. BARKEY **** 

C **** DOUGLAS AIRCRAFT COMPANY , JANUARY 27,1987 **** 

C 

C SUBROUTINES NEEDED: 

C MATRIX-SOLVES SYSTEM OF ALGEBRAIC EQUATIONS GIVEN 

C THE MATRIX OF COEFFICENTS 

C CMSXXX - ALLOWS THE EXECUTION OF CMS COMMANDS INSIDE 

C A FORTRAN PROGRAM 

C ROOT-SOLVES FOR ROOTS OF THIRD OR FOURTH ORDER POLYNOMIALS 

C 

C OUTPUT SENT TO UNIT NUMBER 13, WHICH MUST BE DEFINED AS 

C "PSTRESS OUTPUT" BY A FILEDEF COMMAND. USER WILL BE PROMPTED 

C FOR THE PRINTING OF OUTPUT. 

C 

C W1 AND W2 ARE DEFLECTIONS OF SKIN AND STRINGER FLANGE, 

C RESPECTIVELY. U1 AND U2 ARE AXIAL DISPLACEMENTS OF 

C SKIN AND STRINGER FLANGE, RESPECTIVELY. EP1 AND EP2 

C ARE STRAINS IN SKIN AND STRINGER FLANGE, RESPECTIVELY. 

C TAUC IS SHEAR STRESS IN ADHESIVE. SIGC IS PEEL STRESS IN 

C ADHESIVE. 

C 

C FOR CASES WITH NO AXIAL FORCE: 

C 

C Wl=( Cl*EXP( MU1*X ) +C2*EXP( -MU1*X) +C3*EXP(MU2*X ) +C4*EXP ( -MU2*X) 

C +C5*EXP(MU3*X)+C6*EXP( -MU3*X)+C7*X**3+C8*X**2+C9*X+C10 ) /D1 

C W2=(U21*(C1*EXP(MU1*X)+C2*EXP( -MU1*X) ) +U22* ( C3*EXP(MU2*X) 

C +C4*EXP( -MU2*X) ) +U23* ( C5*EXP(MU3*X)+C6*EXP{ -MU3*X) )+A/B* 

C (C7*X**3+C8*X**2+C9*X+C10) ) /D2 

C U1=(U51* ( C1*EXP ( MU1*X) -C2*EXP( -MU1*X) ) /MU1+U52* ( C3*EXP(MU2*X) - 

C C4 *EXP ( -MU2 *X ) )/MU2+U53*(C5*EXP(MU3*X)-C6*EXP(-MU3*X) ) /MU3- 

C 3* ( TT1+TT2*A/B ) *C7*X**2/ ( F+G) +C11*X+C13 ) / ( E1*T1 ) 

C U2= ( -U51* ( C1*EXP ( MU1*X) ~C2*EXP ( -MU1*X) ) /MU1-U52* ( C3*EXP ( MU2*X ) 

C C4*EXP(-MU2*X) ) /MU2-U53* ( C5*EXP(MU3*X) -C6*EXP( -MU3*X) ) /MU3 

C +3* ( TT1+TT2*A/B ) *C7*X* *2/ ( F+G) +C12*X+C14 ) / ( E2*T2 ) 

C EP1= ( U51* ( C1*EXP ( MU1*X ) +C2*EXP ( -MU1*X ) )+U52* ( C3*EXP(MU2*X) 

C +C4*EXP( -MU2*X) )+U53*(C5*EXP(MU3*X)+C6*EXP(-MU3*X) )- 

C 6* ( TT1+TT2 * A/B ) *C7 *X/ ( F+G ) +C11 ) / ( El+Tl ) 

C EP1=(-U51*(C1*EXP(MU1*X)+C2*EXP(-MU1*X) ) -U52* ( C3*EXP(MU2*X) 

C +C4*EXP(-MU2*X) ) -U53* { C5*EXP(MU3*X) +C6*EXP( -MU3*X) ) 

C +6* ( TT1+TT2*A/B ) *C7*X/ ( F+G) +C12 ) / ( E2+T2 ) 

C 

C NOTE : 

C FOR CASES WITH AXIAL FORCE, SUBSTITUTE EXP(MU4*X) AND 

C EXP(-MU4*X) TERMS FOR X**3 AND X**2 TERMS. 

C 
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c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 


50 

c 

c 


c 


SIG=EC*(Wl-W2)/TC 

TAU=GC* ( 2* (U1-U2 )+Dl*TTl*DWl/DX+D2*TT2*DW2/DX) / ( 2*TC) 

NOTE : C12=( 2* ( TT1+TT2*A/B ) *C8+F*C11 ) /G 

C14=( 6* ( TT1+TT2*A/B ) *C7/ ( E* ( F+G) )+(TTl+TT2*A/B) *C9+ 
F*C13)/G 


CONSTANT 
A 
B 
E 
F 
G 
Cl 
UIJ 
MU I 
T1 , T2 
Dl , D2 
El , E2 
EC , GC 
TC 


PROGRAM EQUIVALENT 
AA 
BB 
EE 
FF 
GG 
E ( I ) 

U(I,J) 

MU ( I ) 

THICKNESS OF SKIN AND STRINGER FLANGE 
RIGIDITY OF SKIN AND STRINGER FLANGE 
MODULUS OF SKIN AND STRINGER FLANGE 
EXTENSIONAL AND SHEAR MODULUS OF ADHESIVE 
THICKNESS OF ADHESIVE 


CHARACTER ANSWER 
CHARACTER* 4 DNAME ( 2 ) 

REAL* 16 D1,D2,E1,E2,T1,T2,TC,GC,EC,L,P1,M1,F1,P2,M2,F2 

REAL* 16 AA , BB , CC , DD , EE , FF , GG , TTl , TT2 , P , Q , R , S , RSIG1 , RSIG2 

REAL* 1 6 RT , X , RU1 , RU2 , RW1 , RW2 , RSIG , RTAU , FLAG 

COMPLEX* 3 2 LAM,U,MU , A,B ,D,E ,F,G,W1 ,W2 ,U1 ,U2 

COMPLEX* 3 2 SIG,SIG1 , SIG2 , TAU 

INTEGER J , N , I , M , COUNT 

DIMENSION LAM( 4 ) ,U( 5 , 4 ) ,MU( 4 ) , A( 12 , 12 ) ,B ( 12 ) ,E( 14 ) ,X( 20 ) 
DIMENSION D(12,0:12) ,W1(20) ,W2(20) ,U1(20) ,U2(20) ,SIG(20) 
DIMENSION SIG1 ( 20 ) ,5102(20) ,TAU(20) ,RU1(20) ,RU2(20) ,RW1(20) 
DIMENSION RW2 (20) / RSIG(20) ,RSIG1( 20) ,RSIG2(20) ,RTAU(20) 

CALL CMSXXX( ' CLRSCRN * ,IRCODE) 

WRITE ( 6 , FMT= ' (//)') 

WRITE (6, 50) 

WRITE (13,50) 

FORMAT (T2 4 , 1 P EEL STRESS’,//) 

**** VALUES OF X FOR PLOTTING: ALL THESE VALUES **** 

**** WILL BE MULTIPLIED BY FLANGE WIDTH/1.2 **** 

DATA (X( I) ,1=1,7 )/0, .01875, .0375, .075, .1125, .15, .225/ 

DATA (X(I) ,1=8, 14)/. 3, . 42 , . 54 , . 66 , . 78 , . 9 , .975/ 

DATA (X( I) ,1=15, 20) /l. 05, 1.0875, 1.125, 1.1625, 1.18125, 1.2/ 
**** INPUT PROPERTIES OF BEAMS AND ADHESIVE **** 

WRITE (6, 100) 

READ* , El , Dl , T1 
WRITE ( 13 , 150 )E1 ,D1 ,T1 
WRITE( 6 , 200 ) 

READ* , E2 , D2 , T2 
WRITE (13, 250 )E2,D2,T2 
WRITE ( 6,300) 


- 27 - 



100 

150 

200 

250 

300 

350 

400 

450 

500 

550 

600 

650 


C 


1500 


C 


800 

C 




SKIN: El, Dl, Tl' ) 

' , 1P1E10 . 3 , OPlFlO . 1 , IX, 0P1F8 . 4 ) 
FLANGE: E2 , D2 , T2 ’ ) 

PROPERTIES OF FLANGE :', 1P1E10 . 3 , OPlFlO . 1 , IX , 0P1F8 . 4 ) 
INPUT ADHESIVE PROPERTIES: EC, GC, TC ’ ) 

ADHESIVE PROPERTIES: ’ , 1P1E10 . 3 , IX , 1P1E10 . 3 , 0P1F8 . 4 ) 

INPUT FLANGE WIDTH: L') 

INPUT SHEAR, MOMENT, AND FORCE APPLIED TO SKIN’ ) 
i SHEAR MOMENT, AND FORCE ON SKIN: , 3F9.2) 

INPUT SHEAR, MOMENT, AND FORCE APPLIED TO FLANGE ) 
SHEAR, MOMENT, AND FORCE ON FLANGE: ', 3F9 . 2 ) 


* * * * 


READ* , EC , GC , TC 
WRITE (13,350) EC , GC , TC 
WRITE ( 6,400) 

READ* , L 

WRITE ( 1 3 , 4 5 0 ) L 
FORMAT (’ INPUT PROPERTIES OF 
FORMAT (’ PROPERTIES OF SKIN: 

FORMAT ( ' INPUT PROPERTIES OF 
FORMAT ( 

FORMAT ( 

FORMAT ( 

FORMAT ( 

FORMAT ( 

FORMAT ( 

FORMAT ( / , 

FORMAT ( ' 

FORMAT ( ' 

M=1 
N=12 
COUNT=0 

**** DEFINE CONSTANTS 
TT1= ( Tl+TC ) /Dl 
TT2=(T2+TC)/D2 
AA=EC/(TC*D1) 

BB=EC/(TC*D2) 

CC=GC*D1*TT1/ ( 4 . *TC ) 

DD=GC*D2*TT2/ ( 4 . *TC ) 

EE=GC/(2.*TC) 

FF=2./(E1*T1) 

GG=2./(E2*T2) 

COUNT=COUNT+l 
U( 1 , 4 ) =0 
U( 2 , 4 ) =0 
U( 3 , 4 ) =0 
U( 4 , 4 )=0 
U( 5 , 4 ) =0 

**** INPUT LOADING ON SKIN AND STRINGER FLANGE 
WRITE (6,500) 

READ*, PI, Ml, FI 
WRITE (13, 800) COUNT 
WRITE (13,550)P1 ,M1 , Fl 
WRITE ( 6,600) 

READ* , P2 ,M2 , F2 
WRITE ( 13,650)P2,M2,F2 

EQOATICEJ 

& R=i^AA+BB+Fl*F2/Dl/D2 ) *EE* ( FF+GG) + ( CC+DD ) * ( AA*TT2+BB*TT1 ) 


**** 
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& +BB*F1/D1+AA*F2/D2 ) 

S=EE* ( FF+GG) * ( AA*F2/D2+BB*F1/D1 ) 

IF (FI .NE. 0 .OR. F2 .NE. 0) THEN 
FLAG=0 
ELSE 

FLAG=1 
END IF 

c **** SOLVE CHARACTERISTIC EQUATION **** 

CALL ROOT ( P , Q , R , S , LAM ) 
c **** CALCULATE EIGENVECTORS **** 

DO 1000 1=1, 4-FLAG 
U(1,I)=1. 

U ( 2 , I ) = ( AA* ( CC+DD ) -F1/D1*DD*LAM( I ) +LAM ( I ) *LAM( I ) *DD ) / 

& ( BB* ( CC+DD ) -F2/D2 *CC*LAM( I ) +LAM ( I ) *LAM( I ) *CC ) 

U ( 3 , I ) =LAM( I ) 

U( 4 , I ) =LAM( I ) *U( 2 , I ) 

U( 5 , I ) =EE*LAM( I ) * (TT1+TT2*U( 2 ,!))/( LAM( I ) -EE* ( FF+GG) ) 

MU ( I ) =CQSQRT ( LAM ( I ) ) 

1000 CONTINUE 

IF ( D1 .EQ. D2) THEN 
U( 2 , 1 ) =-l . 

U( 2 , 2 ) =-l . 

END IF 

C **** DEFINE MATRIX OF COEFFICIENTS FOR BOUNDARY CONDITIONS **** 

C **** NOTE: A*E=B **** 

DO 1050 1=1,12 
DO 1025 J=1 , 12 
A( I , J)=( 0,0) 

1025 CONTINUE 

1050 CONTINUE 

DO 1100 1=1, 4-FLAG 
A ( 5 , 2 * I - 1 ) =MU ( I ) 

A( 5 , 2*1 ) =-MU { I ) 

A(6, 2*1-1) =MU ( I ) *U( 2 , I ) 

A( 6 , 2*1 ) =-MU( I ) *U( 2 , I ) 

A(3, 2*1-1) =U (5,1)* CQEXP ( MU ( I ) *L ) 

A ( 3 , 2 * I ) =U (5,1)* CQEXP ( -MU ( I ) *L ) 

A( 11, 2*1-1 )=-A( 3, 2*1-1) 

A( 11 , 2*1 ) =-A( 3,2*1) 

A(7, 2*1-1) =LAM( I ) * CQEXP ( MU ( I ) *L) 

A( 7 , 2*1 )=LAM( I ) *CQEXP( -MU( I ) *L) 

A( 8 , 2*1-1 ) =A( 7, 2*1-1) *U (2, I) 

A( 8 , 2*1 ) =A( 7, 2*1) *U (2, I) 

A ( 4 , 2 * I - 1 ) =MU ( I ) *LAM ( I ) 

A( 4 , 2*1 ) =-MU ( I ) *LAM( I ) 

A ( 1 0 , 2 * I - 1 ) =U ( 2 , 1 ) 

A( 10, 2*1 )=U( 2,1 ) 

A(l, 2*1-1 )=(CC*( (TT1+F1/D1/CC+TT2*U(2,I) ) *MU( I ) +( FF+GG) * 

& U( 5 , I ) /MU ( I ) ) -MU ( I ) *LAM( I ) ) *CQEXP ( MU ( I ) *L ) 

A( 1 , 2*1 ) =-A( 1 , 2*1 -1 ) *CQEXP ( -2*MU ( I ) *L ) 

A(2,2*I-1)=( DD* ( ( TT1+ ( TT2+F2/D2/DD ) *U ( 2 , I ) ) *MU( I )+( FF+GG) * 
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1100 


c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 


& U(5,I)/MU(I) )-MU(I)*LAM(I)*U(2,I) ) *CQEXP(MU( I ) *L) 

A( 2 , 2*1 ) =— A ( 2 , 2*1-1 ) *CQEXP( -2*MU( I ) *L) 
A(12,2*I-1)=U(5,I)/MU(I) 

A( 12 , 2*1 )=-U (5,1) /MU ( I ) 

A( 9 , 2*1-1 )=-U (5,1) /MU ( I ) 

A(9,2*I)=U(5,I)/MU(I) 

CONTINUE 
A{ 5 , 9 ) =1 . 

A( 6 , 9 ) =AA/BB 

A(3,7)=A(3,7)*( 1-FLAG) -6 . * ( TT1+TT2 * AA/BB ) *L/ ( FF+GG ) *FLAG 
A( 3 , 11 ) =1 . 

A( 11 , 7 )=A( 11 , 7 ) * ( 1-FLAG) +6 . * ( TT1+TT2*AA/BB ) *L/ ( FF+GG) *FLAG 
A ( 1 1 , 8 ) = A ( 1 1 , 8 ) * ( 1-FLAG ) +2 . * ( TT1+TT2* AA/BB ) /GG*FLAG 
A ( 1 1 , 1 1 ) =FF / GG 

A( 7 , 7 ) =A( 7 , 7 ) * ( 1-FLAG ) +6 . *L*FLAG 

A( 7 , 8 ) =A( 7 , 8 ) * ( 1 -FLAG ) +2 . *FLAG 

A( 8 , 7 ) =A( 8 , 7 ) * ( 1 -FLAG ) +6 . * AA*L/BB*FLAG 

A( 8 , 8 ) =A( 8 , 8 ) * ( 1-FLAG ) +2 . *AA/BB*FLAG 

A ( 4 , 7 ) =A ( 4 , 7 ) * ( 1 -FLAG ) + 6 . * FLAG 

A( 10 , 10 ) =AA/BB 

A ( 1 , 7 ) =A ( 1 , 7 ) * ( 1-FLAG) -6 . * ( 1+CC* (TT1+TT2* AA/BB) / ( EE* ( FF+GG) ) ) 

& *FLAG 

A(1 9)=F1/D1 

A(2^7)-A(2,7)*( 1 -FLAG ) - 6 . * ( AA/BB+DD* ( TT1+TT2 * AA/BB ) / ( EE* 

& (FF+GG) ) )*FL AG 

A( 2 , 9 ) =F2/D2* AA/BB 
A( 12 , 12 ) =1 . 

A ( 9 , 7 ) =A ( 9 , 7 ) * ( 1 -FLAG ) +6 . * ( TT1+TT2 * AA/BB ) / ( EE*GG* ( FF+GG ) ) *FLAG 
A ( 9 , 9 ) = ( TT1+TT2 * AA/BB ) /GG 
A( 9 , 12 ) =FF/ GG 

**** BOUNDARY CONDITIONS **** 


**** x = 0 
B ( 5 ) =0 

**** x =0 

B(6)=0 

* ** * x=L 
B ( 3 ) =F1 
**** x=L 
B(11)=F2 
** ** x=L 
B(7)=M1 
**** x=L 
B(8)=M2 
**** x = 0 

B ( 4 ) =0 
**** X = 0 
B ( 10 ) =0 

* * * * X— L 
B ( 1 ) =Pl 
** * * x=L 
B ( 2 ) =P2 


DW1/DX=0 **** 

DW2/DX=0 **** 
E1*T1* DU1 / D X=F1 **** 
E2*T2*DU2/DX=F2 **** 
D1*D**2W1/DX**2=M1 **** 
D2*D**2W2/DX**2=M2 **** 
D**3W1/DX**3=0 **** 

W2=0 **** 

SHEAR CONDITION **** 
SHEAR CONDITION **** 
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c 


Ul— 0 **** 


c 

c 

c 

c 

2230 

2250 

C 

C 


2095 

2100 

2200 

2300 

C 

2805 

2820 

2830 


**** X— o 
B( 12 ) =0 
**** X=0 U2=0 **** 

B( 9 )=0 

**** SOLVE FOR E **** 

CALL MATRIX ( A , B , E , M , N ) 

**** PLUG E INTO A*E=B , TO CONFIRM E **** 

DO 2250 1=1 , 12 
DO 2230 J=l,12 

D ( I , J ) =A ( I , J ) *E ( J ) 

CONTINUE 
CONTINUE 
E( 13 )=E( 12 ) 

**** SOLVE FOR C12 AND C14 **** 

E( 12 )={ 2 . * ( TT1+TT2 * AA/BB ) *E ( 8 ) *FLAG+FF*E ( 1 1 ) ) /GG 
E( 14 ) = ( 6 . * ( TT1+TT2 * AA/BB ) *E { 7 ) *FLAG/ ( EE* ( FF+GG ) ) 

& + ( TT1+TT2* AA/BB ) *E ( 9 ) +FF*E ( 13 ) ) /GG 

**** PRINT OUTPUT TO FILES AND SCREEN (IF DESIRED) **** 
CALL CMSXXX( 'CLRSCRN ' ,IRCODE) 

WRITE (13, 2100) ( I , LAM ( I ) , I , MU ( I ) ,1=1,4) 

WRITE ( 13,2200) ( ( I , J ,U( I , J) , J=1 , 4 ) , 1=1 , 5 ) 

WRITE (13,2300) (I,E(I), 1=1,14) 

ANSWER= ' N ’ 

WRITE ( 6 , FMT= 1 (///,' ' DO YOU WANT TO SEE EIGENVALUES, 

& ' ’EIGENVECTORS, AND COEFFICIENTS ( Y/N) ?’’,/ , 

& ’ ’ "NO" IS DEFAULT. ’ ’ ) ’ ) 

READ (5,5300, END= 2095) ANSWER 
CALL CMSXXX ( 'CLRSCRN ’,IRCODE) 

IF (ANSWER .NE. ’Y’) GOTO 2805 
WRITE (6,2100) (I,LAM(I) ,I,MU(I) ,1=1,4) 

WRITE (6, 2200) ( (I, J,U(I, J) , J=l,4) ,1=1,5) 

WRITE (6,2300) (I,E(I), 1=1,14) 

FORMATt/,’ EIGENVALUES AND THEIR SQUARE ROOTS:', 

& ( 4 ( / , ' LAM( ' ,S,I1, ’ ) ’ , IX, 1P1E15 . 7 , SP , 1P1E15 . 7 , ' i’ , 

& /, ’ MU ( 1 , S , I 1 , ' ) ’ , 2X, 1P1E15 . 7 ,SP, 1P1E15 . 7 , * i * ,/) ) ) 

FORMAT ( ' EIGENVECTORS : ’ , 5 ( / , 

& 4 ( ' U( ’ ,S,I1, ’ , ' ,S,I1, ’ ) ’ ,1X,1P1E15.7,SP,1P1E15.7, ’i’ ,/) ) ) 
FORMAT ( ' COEFFICIENTS ARE: ' ,/, 

Sc 14 ( 5 C( ' ,S,I2, ' ) ' , 2X, 1P1E15 . 7 , SP , 1P1E15 . 7 , 'i' ,/) ) 

**** SUM ROWS OF D AND THEN COMPARE WITH B **** 

DO 2830 1=1,12 
D( I , 0 )=B( I ) 

DO 2820 J=1 , 12 

D ( I , 0 ) =D (1,0) -D (I, J) 

CONTINUE 
CONTINUE 
DO 2900 1=1,12 

IF ( CQABS ( D ( I , 0 ) ) .GT. l.Q-20) THEN 
D ( I , 0 ) =-D (1,0) +B ( I ) 

WRITE (6, 3000)1, B(I) ,D(I,0) 
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WRITE( 13,3000) I ,B ( I ) ,D( I , 0 ) 

END IF 
2900 CONTINUE 

3000 FORMAT (' BOUNDARY CONDITION \S,I2, 

& ' NOT SATISFIED BY COEFFI CENTS .',/,' CONDITION EQUALS:', 

& 1P1E15^7,SP,1P1E15.7, 1 i * , / , ' SHOULD EQUAL: ' ,S, 

& 1P1E15 . 7 SP 1P1E15 . 7 1 i ' ) 

c **** CALCULATE DISPLACEMENTS AND STRESSES FOR VARIOUS **** 

C **** VALUES OF X, AND PRINT THEM **** 

DO 3100 1=1,20 

X( I ) =X( I ) *L/ 1 . 2 
3100 CONTINUE 

DO 4000 1=1,20 
W1(I)=0 
W2 ( I ) =0 
SIGl ( I ) =0 
SIG2 ( I ) =0 
U1(I)=0 
U2 ( I ) =0 
TAU ( I ) =0 

DO 3200 J=l, 4-FLAG 

F=E(2*J-1)*CQEXP(MU( J)*X(I) ) 

G=E ( 2* J ) *CQEXP( -MU ( J ) *X ( I ) ) 

W1 ( I ) =W1 ( I ) +F+G 

W2 ( I ) =W2 ( I ) + ( F+G ) *U ( 2 , J ) 

SIGl ( I ) =SIG1 ( I ) + ( F+G ) *U( 5 , J ) 

SIG2 ( I )=SIG2 ( I ) - ( F+G) *U( 5 , J) 
U1(I)=U1(I)+(F-G)*U(5,J)/MU( J) 

U2 ( I ) =U2 ( I ) - ( F-G ) *U ( 5 , J ) /MU ( J ) 

TAU ( I ) =TAU ( I ) +EE* { ( FF+GG ) *U( 5 , J ) /MU( J ) 

Sc +( TT1 +TT2* AA/BB )*MU(J) )*(F-G) 

3200 CONTINUE 

RW1 ( I ) =REAL ( W1 ( I ) +E( 7 ) *FLAG*X( I ) **3 . +E( 8 ) *FLAG*X( I ) *X{ I ) 

Sc +E(9)*X(I)+E(10) )/Dl 

RW2 ( I )=REAL( W2 ( I‘)+AA/BB* ( E( 7 ) *FLAG*X( I ) **3 .+E ( 8 ) *FLAG 
& *X(I)*X(I)+E(9)*X(I)+E(10) ) )/D2 

RSIGl ( I ) =REAL ( SIGl ( I ) -6 . * ( TT1+TT2* AA/BB ) *E ( 7 ) *FLAG*X ( I ) / 

& ( FF+GG ) +E ( 11 ) )/Tl 

RSIG2 ( I ) =REAL ( SIG2 ( I ) +6 . * ( TT1+TT2 * AA/BB ) *E ( 7 ) *FLAG*X ( I ) / 

Sc ( FF+GG) +E ( 12 ) )/T2 

RU1 ( I ) =REAL ( U1 ( I ) -3 . * ( TT1+TT2* AA/BB ) *E ( 7 ) *FLAG*X { I ) *X( I ) / 
Sc ( FF+GG ) +E ( 1 1 ) *X ( I ) +E ( 1 3 ) ) / ( El *T1 ) 

RU2 ( I ) =REAL ( U2 ( I ) +3 . * ( TT1+TT2 * AA/BB ) *E ( 7 ) *FLAG*X ( I ) *X{ I ) / 
Sc ( FF+GG ) +E( 12 ) *X( I ) +E( 14 ) )/(E2*T2) 

RSIG ( I )=( RW1 ( I ) -RW2 ( I ) ) *EC/TC 

RTAU ( I ) =REAL ( TAU ( I ) +EE* ( FF* ( - 3 . * ( TT1+TT2 * AA/BB ) *E { 7 ) *FLAG 
Sc *X ( I ) *X ( I ) / ( FF+GG ) +E ( 1 1 ) *X ( I ) +E (13) )-GG*(3.* 

Sc ( TT1+TT2 * AA/BB ) *E ( 7 ) *FLAG*X ( I ) *X( I ) / ( FF+GG ) +E ( 1 2 ) 

Sc *X( I )+E( 14 ) ) + ( TT1+TT2 * AA/BB ) * ( 3 . *E ( 7 ) *FLAG*X ( I ) 

Sc *X(I)+2.*E(8)*FLAG*X(I)+E(9) ) ) ) 

4000 CONTINUE 
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4200 

4300 


4500 


4700 

4800 


5000 

5100 


5200 

5300 

6000 


6100 

6200 

6300 

6400 

6500 

7000 


C 


WRITE (6,4200) 

WRITE (13,4200) 

WRITE (6,4300) 

WRITE (13,4300) 

FORMAT (/, 3 4H DEFLECTIONS) 

FORMAT ('5 1H X W1 W2 U1 U2) 

WRITE ( 6,4500) (X(I) ,RW1 ( I ) ,RW2 ( I ) ,RU1 ( I ) ,RU2(I) , 1=1,20) 

WRITE (13,4500) ( X ( I ) , RW1 ( I ) , RW2 ( I ) , RU1 ( I ) , RU2 ( I ) , 1=1,20) 

FORMAT (20(/,0PlF11.4,4(lPlE11.3) ) ) 

WRITE( 6,4700) 

WRITE (13,4700) 

WRITE (6,4800) 

WRITE( 13,4800) 

FORMAT ( / , 3 2H STRESSES ) 

FORMAT ( 52H X SIGC TAUC SIG1 SIG2 ) 

WRITE ( 6,5000) ( X ( I ) ,RSIG( I ) ,RTAU( I ) ,RSIG1 ( I ) , RSIG2 ( I ) , 

& 1=1,20) 

WRITE (13,5000) (X(I) , RSIG ( I ) , RTAU ( I ) ,RSIG1 ( I ) ,RSIG2 ( I ) , 

& 1=1,20) 

FORMAT( 20 ( / , 0P1F11 . 4 , 4 ( 1P1E11 . 3 ) ) ) 

WRITE (6,5200) 

REWIND ( 5 ) 

READ ( 5 , 5300 , END=6000 ) ANSWER 
IF (ANSWER .EQ. 'Y') THEN 

CALL CMSXXX( ' CLRSCRN ',IRCODE) 

GOTO 1500 
ENDIF 

FORMAT (/,' ANOTHER RUN WITH DIFFERENT LOADS ( Y/N )?' ) 

FORMAT (Al) 

WRITE (6, 6200) 

REWIND ( 5 ) 

READ (5,6300, END=7 000) ANSWER 
IF (ANSWER .NE. 'Y') GOTO 7000 
WRITE ( 6,6400) 

REWIND ( 5 ) 

READ ( 5 ,6500, END= 6100) DNAME 
CALL CMSXXX ( ’ EXEC ' , ' DSPR 
& ’A ' , ' ( 

CALL CMSXXX ( 'CP ' , ' SLEEP 

GOTO 7000 

CALL CMSXXX ( ’ EXEC ’ , ' DSPR 
CALL CMSXXX ( 'CP ' , ' SLEEP 

FORMAT ( / , ' SEND RESULTS TO PRINTER ( Y/N )?' ) 

FORMAT (Al) 

FORMAT ( / , ' ENTER PRINTER DESTINATION: < ENTER > FOR DEFAULT') 

FORMAT ( 2 A4) 

CALL CMSXXX ( 'CLRSCRN ',IRCODE) 

STOP 

END 

SUBROUTINE ROOT( P , Q , R , S , LAM) 

**** SOLVES FOR ROOTS OF CHARACTERISTIC EQUATION (LAMBDAS) ie ** i 


' , ' PSTRESS ' , ' OUTPUT ' , 

’ , DNAME, IRA) 

’ , IRA) 

',' PSTRESS ’,' OUTPUT ',IRCC 

' , IRA) 
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REAL* 16 A,B,C, AP , BP , RT, P ,Q , R , S 
COMPLEX* 3 2 L , M , N , AS , BS , LAM 
DIMENSION LAM( 4 ) 

IF (S .EQ. 0) THEN 
A=P 
B=Q 
C=R 
ELSE 

A=(Q-3.*P*P/8. )/2 

B=( (Q-3.*P*P/8. ) * (Q-3 . *P*P/8 . ) -4 . * ( -3 . *P*P*P*P/256 . 

& +P*P*Q/16-R*P/4.+S) )/16. 

C=- ( P*P*P/ 8 . -P*Q/2 . +R) * ( P*P*P/8 . -P*Q/2 . +R) /64 . 

END IF 

AP=B-A*A/ 3 . 

BP= ( 2 . *A*A*A-9 . *A*B+27 . *C ) /27 . 

AS=-BP/2.+CQSQRT(QCMPLX(BP*BP/4.+AP*AP*AP/27. ) ) 

AS=CQABS ( AS ) * * ( 1 . / 3 . ) *CQABS ( AS ) /AS 

BS=-BP/ 2 . -CQSQRT ( QCMPLX ( BP*BP/ 4 . +AP*AP*AP/27 . ) ) 

BS=CQABS(BS)**(l./3. ) *CQABS ( BS ) /BS 

RT=3 . 

L=- . 5 * ( AS+BS ) + ( 0 , 1 . ) *QSQRT ( RT ) /2 . * ( AS-BS ) -A/ 3 . 

M=- . 5* (AS+BS ) - ( 0 , 1 . ) *QSQRT(RT) /2 . * (AS-BS ) -A/3 . 

N=AS+BS-A/ 3 . 

IF (S .EQ. 0) THEN 
LAM ( 1 ) =L 
LAM ( 2 ) =M 
LAM( 3 ) =N 
LAM( 4 ) =0 
ELSE 

L=CQSQRT ( L ) 

M=CQSQRT ( M ) 

N=CQSQRT(N) 

BP=p*p*p/8 . -P*Q/2 . +R 

LAM( 3 ) =QABS ( BP ) /BP* ( -L-N-M) -P/4 . 

LAM ( 2 ) =QABS ( BP ) /BP* ( L-M+N ) -P/ 4 . 

LAM ( 1 ) =QABS ( BP ) /BP* ( -L+M+N ) -P/ 4 . 

LAM ( 4 ) =QABS ( BP ) /BP* ( L+M-N ) -P/4 . 

ENDIF 

RETURN 

END 

SUBROUTINE MATRIX ( D , C , B , M , N ) 

C PURPOSE 

C THIS SUBROUTINE SOLVES A SYSTEM OF LINEAR ALGEBRAIC 

C EQUATIONS BY TRIPLE FACTORIZATION WITH PARITAL PIVOTING 

C AND BACK SUBSTITUTION. 

C 

C TRIPLE MATRIX FACTORIZATION 

C * DECOMPOSE A COMPLEX MATRIX INTO ITS COMPONENT 

C * MATRICES IN THE MANNER A=LDU , WHEREIN 

C * A IS THE GIVEN MATRIX 

C * L IS A LOWER TRIANGULAR MATRIX OF UNIT DIAGONAL 
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C 

* 

D 

IS A DIAGONAL MATRIX 


c 

* 

U 

IS AN UPPER 

TRIANGULAR MATRIX OF UNIT DIAGONAL 

c 

* 

DESCRIPTION 



c 

* 

THE 

ORDER OF DECOMPOSITION 

1 IS IN THE SEQUENCE OF 

c 

* 

COLUMNS. THE LARGEST ABSOLUTE ELEMENT BELONGING 

c 

* 

TO A 

l ROW AS YET 

NOT DECOMPOSED IS USED AS PIVOT. 

c 

* 

THE 

COMPONENT MATRICES SUPERIMPOSE THE ORIGINAL. 

c 

c 

* 

ARG. 

TYPE I /O/S 

DIMS. 

DEFINITION 

c 

* 

A 

COMPLEX* 3 2 

(N,N) 

INPUT AND DECOMPOSED MATRICES 

c 

* 

B 

COMPLEX* 3 2 

(N,M) 

B ON INPUT, X ON OUTPUT 

c 

* 

KOL 

INTEGER 

N 

ROW ORDER OF DECOMPOSITION 

c 

* 

S 

COMPLEX* 3 2 

N 

SCRATCH ARRAY 

c 

* 

N 

INTEGER 


ORDER OF THE SYSTEM 

c 

* 

M 

INTEGER 


NUMBER OF COLUMNS 

c 

c 

* 

CODED 

BY L. CHAHINIAN, 7/16/80 


c * MODIFIED BY D. BARKEY FOR USE IN PEEL STRESS COMPUTATION 

C * PROGRAM "PSTRESS", 2/20/87 

C 

DIMENSION D(N,N) ,C(N,M) ,B(N,M) 

DIMENSION KOL (20) , S ( 2 0 ) , A ( 2 0 , 2 0 ) 

REAL* 16 R 

COMPLEX* 3 2 A,Q,S,B,C,D 
EQUIVALENCE (Q,R) 

DO 200 1=1, N 
DO 100 J=1,N 

A ( I , J ) =D ( I , J ) 

100 CONTINUE 

DO 150 J=1,M 

B ( I , J ) =C ( I , J ) 

150 CONTINUE 

200 CONTINUE 

DO 1000 1=1 ,N 
1000 KOL ( I ) = I 

NM1 = N-l 
DO 4000 K=1,NM1 

q **** DETERMINE LARGEST ABSOLUTE ELEMENT IN THIS COLUMN **** 

KP1 = K+l 
R = CQABS( A(K,K) ) 

L = K 

DO 1200 I=K,N 

I F ( CQABS ( A ( I , K ) ) . LE . R ) GO TO 1200 

R = CQABS ( A( I , K) ) 

L = I 

1200 CONTINUE 

S ( K) = A(L,K) 

IF ( L . EQ . K) GO TO 1241 
J = KOL(K) 

KOL ( K ) = KOL(L) 

KOL ( L ) = J 
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Q O 


C **** INTERCHANGE ROWS K AND L **** 

DO 1220 1=1, N 
Q = A(K,I) 

A( K, I ) = A(L,I) 

1220 A(L,I)= Q 
1241 Q = A(K,K) 

DO 2000 I=KP1,N 
A(I,K)= A(I,K)/Q 
DO 2000 J=KP1,N 

IF ( CQABS { A( K, J ) ) .LT. l.Q-35) A(K, J)=( 0 . , 0 . ) 

2000 A(I,J)= A(I,J)-A(I,K)*A(K,J) 

DO 4000 J=KP1,N 
4000 A(K,J)= A(K,J)/Q 
C 

C **** SOLVE FOR UNKNOWNS, GIVEN TRIPLE FACTORIZATION **** 

C 

C * SOLVE FOR X IN THE MATRIX EXPRESSION AX=B 

C * WHERE THE COLUMNS OF B ARE GIVEN, AND 

C * WHERE A IS GIVEN IN TERMS OF ITS COMPONENT 

C * MATRICES L, D, AND U, SUCH THAT 

C * A = L D U 

C * L IS A LOWER TRIANGULAR MATRIX OF UNIT DIAGONAL 

C * D IS A DIAGONAL MATRIX 

C * U IS AN UPPER TRIANGULAR MATRIX OF UNIT DIAGONAL 

C * DESCRIPTION 

C * THE COMPONENT MATRICES ARE BROUGHT THROUGH ARRAY A. 

C * THEIR ROWS HAVE BEEN INTERCHANGED PER INFORMATION 

C * ARRAY KOL. FOR EACH REQUIRED COLUMN, THE CONTENTS 

C * ARE FIRST TRANSFERRED INTO STORAGE ARRAY S, AND 

C * REARRANGED BACK INTO THE ORIGINAL COLUMN PER KOL. 

C * A FORWARD PASS OF THAT COLUMN ON MATRIX L YIELDS 

C * Y, WHERE LY=B . THIS IS FOLLOWED WITH THE 

C * COMPUTATION OF Z VIA DZ=Y, AND AT LAST UX=Z 

C * PROVIDES X THROUGH A REVERSE PASS ON U. THE 

C * CONTENTS OF ARRAY B HAVE THUS BEEN REPLACED 

C * WITH X. 

C 

C FOR EACH COLUMN OF B 

DO 10000 J=1,M 

C INTERCHANGE ROWS OF B PER KOL 

DO 6100 1=1, N 
6100 S(I)= B ( I , J ) 

DO 6200 1=1, N 
K = KOL ( I ) 

6200 B ( I , J ) = S ( K ) 

**** FOR EACH COLUMN OF B **** 

**** SOLVE LY=B, FOR Y **** 

IM1 = 1 
DO 8000 1=2, N 
DO 7000 L=1 , IM1 
B( I , J ) =B( I , J ) -A( I , L ) *B( L , J ) 
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7000 CONTINUE 
8000 IM1 = I 

C **** SOLVE DZ=Y , FOR Z **** 

DO 9000 1=1, N 
9000 B ( I , J ) = B(I,J)/A(I,I) 

C **** SOLVE UX=Z , FOR X **** 

I = N 

GO TO 9700 

9100 DO 9300 L=IP1 ,N 

9300 B ( I , J ) = B ( I , J ) -A( I , L ) *B ( L , J ) 

9700 IP1 = I 

I = 1-1 

IF ( I . GT . 0 ) GO TO 9100 
10000 CONTINUE 
RETURN 
END 
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